clear
clear matrix
set more off

* ------------------------------------------------------------------------------
* set paths
global root "[include path here]"
global input "$root/Input"
global output "$root/Output"
global intermediate "$root/Intermediate"
global figures "$root/Figures"

* ______________________________________________________________________________
* EVENT STUDY REGRESSION ESTIMATES AND PLOTS

* focus on single family participants who converted once
use "$input/eventStudyPanel_raw.dta", clear
keep if conversions == 0 | (conversions == 1 & wg == "SF")

keep wuse* year month tau conversions pid

local ewin = 60 // event window
set matsize 500

* merge application information into panel
sort pid year month
merge m:1 pid using "$input/applications.dta"
drop if _merge == 2
tab _merge conversions
drop _merge

* generate event time relative to application date
gen apptau = 12*(year - year(appdate)) + month - month(appdate) if conversions == 1
tab tau if apptau == 0

* make sample as small as possible
compress
sort pid

******* event study relative to enrollment date
* create event dummies
* lag dummies
foreach t of numlist `ewin'/1 {
local id = `t'
gen byte D_n`id' = 0
replace D_n`id' = 1 if tau == -`t'
}
* lead dummies
foreach t of numlist 0/`ewin' {
local id = `t'
gen byte D_p`id' = 0
replace D_p`id' = 1 if tau == `t'
}

* create the bumper dummies
replace D_n`ewin' = 1 if tau < -`ewin' & conversions == 1
replace D_p`ewin' = 1 if tau >  `ewin' & conversions == 1

** regressions relative to enrollment date
areg wuse D_n`ewin'-D_n1 D_p1-D_p`ewin' i.year#i.month if conversions == 0 ///
| (conversions == 1 & inrange(tau, -`ewin', `ewin')), a(pid) vce(cl pid)
estimates save "$output/es_tau_ew`ewin'_pFEymFE.ster", replace

******* event study relative to application date
drop D_*
* create event dummies
* lag dummies
foreach t of numlist `ewin'/1 {
local id = `t'
gen byte appD_n`id' = 0
replace appD_n`id' = 1 if apptau == -`t'
}
* lead dummies
foreach t of numlist 0/`ewin' {
local id = `t'
gen byte appD_p`id' = 0
replace appD_p`id' = 1 if apptau == `t'
}

* create the bumper dummies
replace appD_n`ewin' = 1 if apptau < -`ewin' & conversions == 1
replace appD_p`ewin' = 1 if apptau >  `ewin' & conversions == 1

** regressions relative to application date
areg wuse appD_n`ewin'-appD_n1 appD_p1-appD_p`ewin' i.year#i.month ///
if conversions == 0 | (conversions == 1 & inrange(apptau, -`ewin', `ewin')), ///
a(pid) vce(cl pid)
estimates save "$output/es_apptau_ew`ewin'_pFEymFE.ster", replace

******* plot event study relative to application date
clear matrix
set mat 1050

local tau0 "application" // enrollment
local normOp "1000 gal" // "gal/ft{sup:2}"
local m "es_apptau_ew60" 
local mat "pFEymFE"
local xlines "xline(5, lc(gs5) lpattern(-)) xline(0)"

* load estimates from event study regressions
estimates use "$output/`m'_`mat'.ster"
mat `m'_`mat' = e(b)', vecdiag(e(V))'
svmat `m'_`mat', names(matcol) // convert from matrix to variables

* create 95% confidence intervals
gen `m'_`mat'CIlb = `m'_`mat'y1 - (1.96*sqrt(`m'_`mat'r1))
gen `m'_`mat'CIub = `m'_`mat'y1 + (1.96*sqrt(`m'_`mat'r1))

* Only need the coefficients on the lead/lags
keep if _n <= 120 // the number of D_* coefficients
expand 2 in 1 // create an extra row for the omitted category

* create horizontal axis variable: event time
gen eventTime = _n
replace eventTime = eventTime - 61
replace eventTime = eventTime + 1 if eventTime >= 0
replace eventTime = 0 if eventTime == 61

sort eventTime
order eventTime 

* change to zeros the row values for the omitted cateogry
local suffixname "y1 r1 CIlb CIub"

foreach suffix of local suffixname {
  replace `m'_`mat'`suffix' = 0 if eventTime == 0
}

* ---------------------------------------------------------------
* Plots (all methods)
gr tw sc `m'_`mat'CIlb `m'_`mat'CIub `m'_`mat'y1 eventTime, ///
c(l l i) msymbol(i i smcircle) mcolor(gs10 gs10 gs3) ///
xti("{&tau}: event time in relation to `tau0' date (months)") ///
yti("Estimate of {&kappa}{sub:j} [`normOp']") `xlines' ///
legend(off) ylabel(, nogrid) xlabel(-60(12)60) ///
graphregion(color(white)) lcolor(gs10 gs10 gs3)
graph export "$figures/figure_3a.pdf", replace

clear all
clear matrix